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A theoretical approach was developed for an exact numerical description of a pair of ultracold 
atoms interacting via a central potential that are trapped in a three-dimensional optical lattice. 
The coupling of center-of-mass and relative-motion coordinates is explicitly considered using a 
configuration-interaction (exact-diagonalization) technique. Deviations from the harmonic approxi- 
mation are discussed for several heteronuclear alkali-metal atom pairs trapped in a single site of an 
optical lattice. The consequences are discussed for the analysis of a recent experiment [C. Ospelkaus 
et al, Phys. Rev. Lett. 97, 120402 (2006)] in which radio-frequency association was used to create 
diatomic molecules from a fermionic and a bosonic atom and to measure their binding energies close 
to a magnetic Feshbach resonance. 



I. INTRODUCTION 

The physics of ultracold atoms has attracted a 
lot of interest since the experimental observation of 
Bose-Einstein condensation in dilute alkali-metal atom 
gases [l|, Q- Besides the exciting physics at ultracold 
energies by itself, a further important progress was the 
loading of the ultracold gas into an optical lattice formed 
with the aid of standing light waves 0, 0, 0] • The optical 
lattice resembles in some sense the periodicity of a crystal 
potential @, 0j Si • contrast to real solids the lattice 
parameters are, however, easily tunable by a variation 
of the laser intensity (trap depth) or wavelength (lattice 
geometry). Moreover, different atoms possess different 
interaction potentials that can be either attractive or re- 
pulsive. While different kinds of chemical elements, their 
isotopes, or atoms in different electronic or spin states 
cover already quite some range of interaction strengths, 
a further tunability of the atom-atom interactions in op- 
tical lattices can be achieved using magnetic Feshbach 
resonances [1, [l^l ■ Close to the resonance value of the 
magnetic field the interaction varies in a wide range of 
attractive and repulsive values. 

Ultracold atoms deposited in light crystals are ideal 
systems for a realization of the Hubbard model 0] . This 
model takes into account a single band of a static lat- 
tice potential and assumes the interactions to be purely 
local (low tunneling limit). In this case the optical 
lattice is considered as an array with a very small filling 
rate; optimally with one or two atoms per lattice site. 
The experimental study of a bosonic Mott insulator Q 
and a fermionic band insulator Q provided such a sys- 
tem. The fact that the analysis of a single site is sufficient 
for such systems simplifies their theoretical description 
drastically. 

Another interesting aspect is that ultracold atoms can 
bind together to form ultracold molecules. The optical 
lattice can shield the often fragile, very weakly bound 
molecules from destructive three-body collisions. The 
physics of ultracold atom pairs in optical lattices with 
controllable interactions is thus presently an intensively 



investigated research area. Recently the observation of 
confinement-induced molecules, repulsively interacting 
pairs, and real molecules for both homonuclear [T^l 
and heteronuclear atomic species in optical lattices 
has been reported. 

In order to describe the behavior of atoms in an op- 
tical lattice the latter is usually considered as an array 
of harmonic traps. In such an approach some important 
features of the optical lattice can be lost. For example, 
the correct sinusoidal potential exhibits an energy band 
with a spread of transition energies while the harmonic 
potential possesses a discrete equidistant spectrum. Nev- 
ertheless, the experiment of Stoferle et al. showed 
good agreement with a simplified theoretical description 
based on the harmonic approximation. In their analysis, 
Stoferle et al. compared the measured binding energies of 
confinement-induced molecules and real molecules to the 
ones predicted by a simplified theory in which two atoms 
arc trapped in a harmonic potential and interact via a 
(5-function pseudopotential. Within such a model an an- 
alytical solution exists in the case of two identical atoms 
(in the same quantum states) [T3 |. However, another 
experiment that adopted higher resolution spectroscopy 
and considered a heteronuclear system was interpreted 
to clearly indicate a break-down of the harmonic approx- 
imation p^ . 

From the theoretical point of view, the description of 
two atoms in an optical lattice beyond the harmonic 
approximation is very laborious. The anharmonic part 
of the optical lattice potential leads to a coupling of 
center-of-mass and relative motion and requires there- 
fore to solve the full six-dimensional problem. In fact, 
even within the harmonic approximation different trap- 
ping potentials experienced by the two atoms lead to a 
coupling of center-of-mass and relative motion [lH . [l^. 
This situation occurs, e.g., for heteronuclear atom pairs 
or two atoms of the same kind but in different electronic 
states. 

In this work a numerical approach is developed that 
allows in principle to describe two atoms trapped in an 
optical lattice in an exact way, if the interatomic inter- 
action potential is central (isotropic) and can be given 
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in terms of a single potential curve. The fact that this 
latter curve stems usually from a full molecular calcu- 
lation and is only given in numerical form is explicitly 
considered. Extensions to non-central (like dipolar) in- 
teractions are rather straightforward and planned for the 
near future. The anharmonic coupling is treated in a 
configuration-interaction (CI) like fashion also known as 
exact diagonalization. It leads for sufficiently large ex- 
pansion lengths to exact results. Although the approach 
allows to consider multiple-well potentials and thus more 
than a single site of an optical lattice the present 
work focuses on results obtained for two atoms in a sin- 
gle site. Motivated by the experiment reported in [l3| 
the heteronuclear atom pair formed by fermionic '"^K and 
bosonic ®^Rb was used as a generic system in the present 
study. In order to investigate the influence of the atomic 
interaction strength, its value was varied artificially by 
a controlled manipulation of the inner wall of the cor- 
responding potential curve. For the investigation of the 
influence of the mass difference on the results, two other 
heteronuclear pairs, ^Li-^Li (almost equal masses) and 
6Li_i33Qg (very large mass difference) were considered. 
After a systematic investigation of the effects of anhar- 
monicity and coupling of center-of-mass and relative mo- 
tion a comparison to the experimental data as well 
as to a subsequent theoretical analysis [3 performed in 
parallel to the present work is given. 

The paper is organized in the following way. In SccHIl 
the theoretical approach is described. This includes the 
formulation of the problem in Sec. Ill Al as well as a de- 
scription of the trap parameters (Sec. Ill Bj) . the used in- 
teratomic interaction potentials fSec. lII C"]) . and its sys- 
tematic variation in Sec. Ill Dl The section ends with the 
computational details described in Sec. Ill El The results 
presented in Scc. lIIII are first discussed in terms of en- 
ergies (and their differences) for the generic ^^Rb-^^K 
dimcr in Sec. lIII A|) and for other dimers in Sec. / IIIIBI 
This is followed by an analysis of the radial pair densi- 
ties (Sec. lIII C 1|) and the wavefunctions in absolute co- 
ordinates (Sec. lIII C 2[) . This is followed by a comparison 
to the experimental and recent alternative theoretical re- 
sults in Sec. lIII D\ Finally, a conclusion and outlook is 
given in Sec. lIVI All equations and quantities in this pa- 
per are given in atomic units unless otherwise specified. 



II. SYSTEM 

A. Hamiltonian 

The Hamiltonian describing the interaction of two 
atoms with coordinate vectors ri and f2 trapped in a 
three-dimensional optical lattice is given by 

H(fi,f2) = Ti(7"i) + T2(f2) + U(7"i, f2) 

+ Vti.ap,l(?''l)+Vt,.ap,2(r2) (1) 



where Tj is the kinetic energy operator for particle j, 

U is the atom-atom interaction potential, and Vtrapj is 
the trapping potential for particle j. For optical lattices 
Vtrap.j is often (and also in the present work) given by 

Vtrapj = J2 Sin2(fc,c,), (2) 

c—x^y,z 

where V^f is the potential depth which particle j expe- 
riences along direction c, and kc ~ 2tt/Xc is the wave 
vector and Ac the wavelength of the laser creating the 
lattice potential along the (Cartesian) coordinate c. 

A direct solution of the Schrodinger equation with the 
Hamiltonian given in the form of Eq. ([T]) is complicated, 
since U depends in general on all six coordinates describ- 
ing the two-particle system, even if the atom-atom in- 
teraction is central, i.e. U = U(r) with r = \ri — r2\- 
For realistic interatomic interaction potentials (that are 
usually even only known numerically) , there is no separa- 
bility and this leads to very demanding six-dimensional 
integrals. Therefore, it is more convenient to treat the 
two-particle problem in center-of-mass (COM) and rela- 
tive (REL) motion coordinates. If spherical coordinates 
are adopted, a central interaction potential leads to a 
function of the radial coordinate only. 

On the other hand, the formulation of the two-particle 
problem in COM and REL coordinates complicates the 
treatment of the trap potential, since its separability in 
Cartesian coordinates is lost in the COM and REL coor- 
dinate system. However, performing a Taylor expansion 
of the sinusoidal trapping potential ^ around the origin 
simplifies the problem drastically, because the angular 
parts can be analytically solved for in the case of a cen- 
tral interatomic interaction potential. For two identical 
atoms in the same state the use of the harmonic approxi- 
mation for the trapping potential leads even to a problem 
that is completely separable in COM and REL coordi- 
nates p^ . If the true atom-atom interaction is, further- 
more, replaced by a i5-function pseudopotcntial that re- 
produces only asymptotically the two-body zero-energy 
s-wave scattering, the Schrodinger equation possesses an 
analytical solution for both isotropic or anisotropic har- 
monic traps p^ . [20| . Noteworthy, even within the har- 
monic approximation the separability is lost, if the two 
atoms experience different trapping potentials. This is 
the case, if a heteronuclear system or two identical atoms 
in different electronic states are considered. 

After performing the Taylor expansion of the sinusoidal 
trapping potential ^ around the origin, the transfor- 
mation of the Hamiltonian IH) into the COM and REL 
coordinate systems leads to a Hamiltonian of the form 

H(i?,f) =hcOM(i?) + hREL(r)+W(i?,f) (3) 

with 

hcOM(i?) =tkin(i?)+VOL(i?), (4) 

hREL(r) = tkin(r) + Vol (r)+ Vint (r) • (5) 
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It is worth emphasizing that in the present formulation 
only the truly non-separable terms (represented by prod- 
ucts of COM and REL coordinates) are left in the cou- 
pling term W. All separable terms of the optical lat- 
tice (OL) potential are included into the COM and REL 
Hamiltonians hcoM and Iirel respectively. 

In a first step, the eigenstates and -values of the COM 
and REL Hamiltonians are obtained independently of 
each other by means of a numerical solution of the cor- 
responding stationary Schrodinger equations. 

hcoM \ipt) = £i l-iAj) and hREL = I'f'i) ■ (6) 

The wavefunctions ip{R) and (/)(r) are then used to 
form the configuration state functions 7^) ~ 

'ipii,{R) 4>jf^{f). The stationary Schrodinger equation 
with the full Hamiltonian given in Eq. ([3]), 

H|vl/,) , (7) 

is then solved by expanding ^E* as ^{R,r) = 
J2k ^kiR, r). Insertion of this expansion into Eq. ([7]) 
leads (after the usual manipulations) to a matrix eigen- 
value problem that is solved numerically and yields the 
energies £i and eigenvector coefficients Ck- 

B. Trap parameters 

Despite the already mentioned breakdown of the har- 
monic approximation especially for heteronuclear sys- 
tems it is still convenient to introduce the mean 
harmonic-oscillator frequencies and Jlho of a single 
lattice site for the REL and COM motion respectively, 

^ k ^2 ^^^-^;^-^^^ (8) 

17ho = k ^2^^ • (9) 

In Eqs. ([5]) and ^ fi and M denote the reduced mass 
and total mass of the two particles respectively, Hj is 
defined as fii^2 — ^■/f^2,i where mj is the mass of atom 
j, and Vj = lo-aj is the optical lattice depth that is equal 
to the product of the laser intensity Iq (for an isotropic 
geometry Iq = Ix ~ ly = Iz) and the polarizabilities 
aj of atom j. Finally, one has k = = ky = for 
an isotropic geometry of the lattice. This isotropy is in 
fact assumed in Eqs. ^ and For identical particles 
of mass m Eq. ([5]) reduces to the well-known relation 

coho = k^y2Vo/m 

Some parameters of the trap chosen in the present 
study were motivated by the recent experiment reported 
in [13| . Therein a three-dimensional optical lattice gen- 
erated by lasers with wavelength A = Aj. = Ay = A^ of 
1030 nm was used for the trapping of ultracold bosonic 
^^Rb and fcrmionic ^''K atoms. The two different lattice 



depths Vftb = 40 Ef^ and Vnh = 27.5 E^^ were consid- 
ered where the individual recoil energy is defined, e.g., as 
E^'^ = k"^ / (2mRb)- Since the static dipole polarizabilities 
of rubidium and potassium are different, aRb = 324 a. u. 
and ax = 301a. u. [2l| . the two atoms experience dif- 
ferent potentials: Vk = 37.2 £;^'' and = 25.5 for 
AQEf^ and 27.5 £;^^ respectively. The mean harmonic- 
oscillator frequencies ^ are Who(40 E^"") = 27rx 35.7 kHz 
and Who(27.5 Ef^') = 27r x 30 kHz. While most of the re- 
sults of this work are obtained for these frequencies, some 
other values are also considered in order to investigate the 
influence of the trap frequency in more detail. 



C. Interatomic interaction potential 

The interaction between rubidium and potassium 
atoms is modeled using the Born-Oppenheimcr (BO) po- 
tential of the a'^I]+ electronic state describing the in- 
teraction of two spin-polarized atoms. In general the 
atom-atom interaction potentials are only known numer- 
ically. For the short-range part Vsr of the potential in 
between R e [1.588 ao, 18.2 oq] the data of [231 ^'^^ used 
(ao is the Bohr radius). The data points at i? = 17.6 oq 
and R = 16.99998 oq have been omitted, because their 
inclusion results in a non-smooth potential curve. The 
long range part Vlr of the a'^E+ electronic state is 
constructed in a similar way as was done by Zemke et 
al. [2^. Therefore, this long range part is defined as 
Vi^^ir) = D, + Aydisp(r) + AKx(0 for R > 18.6 
where AKiisp(r) = -Cg/r^ - Cg/r^ - Cio/r^" and the 
dispersion coefficients C„ are the values of Derevianko 
and co-workers [24| .l25l except Cq = 4292 ± 19 a. u. which 
was taken from [26|. The exchange interaction is given 
by AVo^{r) = -Cr'^e-^'' with C = 0.00231382, a = 
5.25603, /? = 1.11892 as given in 0]. To merge the short- 
and the long-range parts the short-range part is raised up 
by half of the value iJmerge = Fsr(18.2 ao) - Vlr(18.6 uq). 
According to [2^ the a'^E+ state supports 32 bound 
states and the interaction of the atoms via the a'^E+ 
potential is strong and repulsive. The same amount of 
bound states and the same character of the interaction 
are observed in the present calculation using the poten- 
tial curve constructed the way described above. 

In order to study the influence of different masses, po- 
larizabilities, different interaction potentials, and also to 
check the generality of the conclusions of this work other 
systems are also analyzed. In particular, the heteronu- 
clear ^Li-'^Li and ^Li-^^'^Cs pairs interacting via their re- 
spective a^T,~^ electronic state are considered. The po- 
tential curves for ^Li-''Li and ^Li-^'^^Cs were constructed 
according to [lB| and [l^l , respectively. 



D. Manipulation of the interatomic interaction 

In the limit of zero collision energy the interaction be- 
tween two atoms can be characterized by their s-wave 
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scattering length Osc- The sign of age determines the type 
of interaction (repulsive or attractive) and the absolute 
value determines the interaction strength. Experimen- 
tally it is difficult to accurately measure the scattering 
length. For example, there is no agreement about the 
value of the triplet scattering length for the ®^Rb-^°K sys- 
tem. According to the ongoing discussion [2^, [H, [2^, [s^l 
the value —185(4) oq appears to be the most reliable one. 
A standard way to match the calculated scattering length 
with the experimental value is a smooth shift of the in- 
ner wall of the BO potential as is described in [2^. This 
procedure can also be used for an effective variation of 
the scattering length since a systematic variation of the 
inner wall allows to shift the least bound (lb) state sup- 
ported by the potential curve as is shown in Fig.[T] If the 
least bound state is close to the dissociation threshold or 
moves even into the dissociative continuum, the scatter- 
ing length and thus the interaction between the ultracold 
atoms are strongly influenced. Therefore, a small vari- 
ation of the inner wall of the potential can modify the 
interatomic interaction potential from strongly repulsive 
to strongly attractive. This procedure is adopted in the 
present work in order to investigate the influence of the 
interatomic interaction potential. The scattering length 
is only well-defined for zero-energy scattering and thus 
the underlying concept is in principle not applicable to 
trapped particles with a non- vanishing zero-point energy. 
Therefore, the scattering-length values (for a given inner- 
wall shift) arc determined for the trap-free situation. In 
this case age can be uniquely determined from the analy- 
sis of the shape of the zero-energy scattering wave func- 
tion [U. 

Ignoring the formal problems of defining a scattering 
length within a trap (that will be discussed in some more 
detail below), it is often considered useful to introduce a 
dimension-free interaction parameter ^ that reflects the 
relative magnitude of the interaction strength with re- 
spect to the confinement by the trap. If this confinement 
is approximated within the harmonic approximation, the 
interaction parameter is naturally defined as ^ = asc/oho 
where Oho is the characteristic length of a harmonic po- 
tential given by Oho = 1 / ^/.io^ho • For a heteronuclear 
atom pair tJho is again the mean harmonic frequency de- 
fined in Eq. 

Experimentally, a strong variation of the interaction 
strength can be realized with the help of magnetic Fesh- 
bach resonances (MFR) d, The MFR technique 

was also used to tune the interatomic interaction from 
strongly repulsive to strongly attractive in the already 
mentioned experiments with atoms in the optical lat- 
tices [H, [12, [3 ■ general, the correct theoretical de- 
scription of a MFR requires a multi-channel scattering 
treatment which in the present case would have to in- 
corporate also the optical lattice. In the analysis of the 
experiments described in [l^, [ll] it is, however, as- 
sumed that it is possible to model the MFR in an ef- 
fective two-channel picture [3l|. Within this model it is 
straightforward to relate the applied magnetic field to a 




FIG. 1: (Color online) Sketch (not to scale) of a cut through the 
potential surfaces along the x direction (y = z = 0) for a system of 
two identical atoms if one of them is positioned at the zero of x. The 
upper graph shows the sin^ potential together with the harmonic 
(green) and sextic (red) approximations. While the harmonic and 
sextic potentials support solely bound states, the energy spectrum 
of the sin^ potential is partly discrete (for a sufficiently deep value 
of Viat) and partly continuous. The lower graph shows the range 
of small X values on an enlarged scale. A very tiny variation of 
the inner wall of the interaction potential leads to a relatively large 
shift of the least bound and the first trap-induced states. 



scattering- length value (see Eq. ([^D|) below) . 



E. Computational details 

The eigenfunctions of the Hamiltonians hcoM and 
liREL are obtained by expressing both ip{R) and 4>{f) 
as a linear combination of products of radial i?-spline 
functions times spherical harmonics. The corresponding 
Schrodinger equations are solved numerically using the 
Rayleigh-Ritz-Galerkin approach [s^ which leads to an 
algebraic eigenproblem. 

In general, the lattice leads to a coupling of the an- 
gular momenta. Therefore, the spherical harmonics are 
no eigensolutions of the angular part. Due to the cu- 
bic trap geometry used in the experiment in and 
also for the present calculations, the coupling of different 
spherical harmonics is weak. In fact, the orbitals ip{R) 
and 4'{r) describing the states relevant to this work are 
almost converged, even if only Z = is considered. How- 
ever, the coupling term W in the Hamiltonian ([3]) leads 
to a stronger angular momentum coupling. Good con- 
vergence was found in the CI calculation, if all spherical 
harmonics up to Z = 3 (and thus also —3 < m < +3) 
were included in the calculation of the orbitals i/j^R ) and 
</>(f). 

The required number of B splines and their knot se- 
quence depend strongly on the behavior of the wave func- 
tion {ip{R) or 0(r)) that should be described. In the 
context of ultracold collisions the main interest is put on 
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the energetically low-lying COM orbitals ^{R ) that pos- 
sess a small number of nodes. For the results discussed in 
this work, about 70 B splines were found to be sufficient 
to obtain convergence. Evidently, more complicated or 
highly anisotropic trap geometries (like double or triple 
wells [l3|) require larger expansions. 

The numerical description of the REL orbitals 4'{f) is 
more demanding, if one is interested in the most weakly 
bound states or the low-lying dissociative states. The BO 
curves of alkali-metal atom dimers support often a large 
number of bound states (for example, the ^^Rb-^°K sys- 
tem possesses in the a^S+ state already 32 bound states 
for / = 0). The very long-ranged, weakly bound states 
consist therefore of a highly oscillatory inner part (cover- 
ing the so-called molecular regime and providing the or- 
thogonality to all lower lying bound states) and a rather 
smooth long-range part. Correspondingly, it is practical 
to use two different knot sequences for the B splines. In 
the present case convergence was found if 200 B splines 
expanded on a linear knot sequence covering the interval 
< r < 20 ao are used together with about 70 B splines 
for the remaining r range. The latter 70 B splines are 
expanded on a knot sequence in which the separation 
between the knot points increases in a geometric fashion. 

Converged CI calculations were found, if they com- 
prised configurations built from about 120 REL and 60 
COM orbitals. After taking symmetry into account this 
amounts to about 1060 configurations forming the CI ex- 
pansion for the states of interest in this work. 



III. RESULTS 
A. Energy spectrum of the ^^^Rb-^^K system 

The description of an optical lattice beyond the har- 
monic approximation is in the present work achieved by 
extending the Taylor expansion of the sin^ potential be- 
yond the harmonic (1st order and thus quadratic) term. 
In principle, one should seek for convergence with respect 
to the expansion length, but there are some practical rea- 
sons why a simple convergence study as a function of the 
expansion length causes problems. First of all, even-order 
expansions like the 2nd order one which leads to polyno- 
mials with a degree of up to 4 (quartic potential) support 
an infinite number of bound states with negative energy, 
since they tend to — oo for x approaching either -|-oo or 
— cx). However, these bound states with negative ener- 
gies are unphysical, since they do not exist in the case 
of the (original) positive definite sin^ potential. The 3rd 
order expansion that leads to polynomials up to a degree 
of 6 (sextic potential) supports on the other hand (like 
all odd-order expansions) only bound states with positive 
energy values. 

A comparison of this sextic potential with the sin^ po- 
tentials shows that the sextic potential reproduces ex- 
tremely well a single site of the sin^ potential and thus 
of the optical lattice (see Fig. [T]). Therefore, the sextic 



potential is a good choice for the investigation of the ef- 
fects of anharmonicity on the bound states in a single 
site of an optical lattice. Evidently, the sextic poten- 
tial cannot reproduce effects that are due to tunneling 
between neighbor potential wells. Therefore, extended 
(energetically higher lying) bound states in the optical 
lattice that are markedly affected by tunneling are not 
well reproduced by a sextic potential. Noteworthy, even 
in this case a simple convergence study will, however, 
not work. For example, the 5th order and thus next 
odd-order expansion shows a triple-well structure, but 
the two outer wells have a depth and width that differs 
pronouncedly from the correct shape (and the central 
well). This leads to completely wrongly described states 
in these outer wells and may thus show wrong tunneling 
behavior for the states in the middle well, especially in 
the case of resonant tunneling. Since the present study 
concentrates on the anharmonicity effects within a single 
site of an optical lattice, only the sextic potential and, 
for comparison, the harmonic one are considered. Effects 
that are due to tunneling between neighbor wells and 
thus include more than single- well potentials are investi- 
gated in a separate work |17| . 

The potential seen by the two atoms in an optical lat- 
tice contains, of course, in addition to the trap poten- 
tial also the interatomic interaction potential that in the 
present case is described by a Born-Oppenheimer poten- 
tial curve (and appears only in the REL coordinates). As 
is sketched in Fig. [1] the interatomic interaction domi- 
nates the short-range part of the potential and leads in 
the case of alkali- metal atoms to a large number of bound 
molecular states. Since the trap potential is compared 
to the variation of the BO curve almost constant in the 
range of the molecular bound states, especially the lower 
lying of these states will in practice not be influenced by 
the optical lattice. The largest possible effect of the op- 
tical lattice on the molecular bound states is expected 
to occur for the energetically highest lying one, the least 
bound (lb) state. 

Due to the large spatial extension of the trap states 
of typical experimentally realized optical lattices these 
trap-induced states are expected to be only weakly influ- 
enced by the molecular potential. However, an immedi- 
ate consequence of the existence of the molecular bound 
states below the trap-induced ones is the nodal structure 
at short distances that is imprinted on the wavefunctions 
and leads to the required orthogonahty of the eigenstates. 
The energetically lowest lying and thus first trap-induced 
(Iti) state possesses thus exactly one more node than the 
lb state. In the experiments most closely related to the 
present work [13, the transition energy between the lb 
and the Iti state has been measured by either rf dissoci- 
ation or association, respectively. This transition energy 
was called binding energy, but it should be kept in mind 
that its definition does not coincide with the standard 
definition of a molecular binding energy which is given 
by the energy difference between a molecular bound state 
and the (lowest) dissociation limit. In the present case 
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the existence of the optical lattice leads to a discretization 
of the dissociation continuum and thus to an additional 
energy shift due to the zero-point energy of the trap. 

If the coupling W of REL and COM coordinates is 
ignored, the energies Ei^ and Em of the least bound 
and the 1st trap-induced state, respectively, are obtained 
from the eigenvalues of Eq. ([6|) as 
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where n specifies the expansion length describing the op- 
tical lattice: n = 2 for a harmonic and n = 6 for a sextic 
trap. In accordance with the underlying assumption of 
an ultracold gas, the system is assumed to be in its low- 
est state with respect to translational motion, i. e. in the 
COM ground state with energy ei. The corresponding 
wavefunctions are given by <i>ib(-R, r) ~ V'i(-R)(/'ib(^ and 

After the inclusion of the coupling of REL and COM 
motion the wavefunctions $ are no longer eigenstates of 
the Hamiltonian, but are used as a basis for expanding 
the full wavefunctions ^' that are obtained together with 
their energies £ by solving the Schrodinger Eq. ([7]) . The 
state 4* with a dominant contribution from $ib ('I'lti) 
is then identified as least bound (Ist-trap-induced) state 
with energy £^^^ (S^), where n stands again for the 
order of the Taylor expansion of the optical-lattice po- 
tential. 

In Fig. [2] the energies of the least bound and the 1st 
trap-induced state are shown for ^^Rb-*°K as a function 
of the trap-free scattering length (see Sec. IIIPp for dif- 
ferent levels of approximation ranging from the separable 
harmonic one to the fully coupled sextic solution. The 
trapping parameters were chosen in accordance with the 
corresponding experiment p^ . Clearly, the energies for 
the different approximations differ most for large posi- 
tive scattering lengths and thus in the case of a strong 
repulsive interaction between the atoms. 

A comparison of the results for the least bound and 
the first trap-induced states reveals that the energy of 
the former is almost unaffected by the anharmonicity of 
the trap and COM-REL coupling. Although this state 
is at least in the strongly repulsive part of the spectrum 
long ranged, it "feels" the anharmonic form of the trap- 
ping potential very weakly. This state remains thus suf- 
ficiently deeply localized in the trap potential and the 
harmonic approximation works still reasonably well; even 
for the rather strong repulsive interaction expressed by 
the scattering length Osc ~ 6500 cq. The energy change 
of the first trap-induced state due to the anharmonicity 
and REL-COM coupling is on the other hand much more 
pronounced. This change is thus predominantly defining 
the modification of the binding energy due to the trap. 

If the results for different levels of approximation are 
compared with each other, a clear ordering is visible. In- 
dependent of the scattering length and thus the interac- 
tion strength as well as its type (repulsive or attractive) 
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FIG. 2: (Color online) Energies of the 1st trap-induced (main 

graph) and least bound (inset) state of ®''Rb-*''K dimers in a single 

site of an optical lattice (A = 1030 nm) for the potential depths 

(a) Vab = 27.5 E^*^ and (b) Vftb = 40.0 E^'' at different levels 

of approximation as a function of the trap-free scattering length 

(see Sec. Ill Dl for details). The energies obtained with a full CI 

calculation for a sextic potential (£(®) , green chain) and a harmonic 

one (f''^', red dashes) are compared to the corresponding sextic 

{E^^\ blue dots) and harmonic {E^'^\ black solid) energies that 

are obtained, if the coupling between COM and REL motion is 
3 

neglected. (Note, the sum -(ti^ho-l-S^ho) corresponds to E^'^^{aBc = 
0) = 100.65 kHz.) 



the uncoupled harmonic energy i?*-^^ is lowered, if the 
COM-REL couphng is included (f (^'). Note, this cou- 
pling exists even within the harmonic approximation for 
a heteronuclear diatomic molecule like RbK, since the 
two atoms possess different masses and polarizabilities 
and experience therefore different trap potentials. As a 
consequence, COM and REL motions do not separate. 
Only for diatomic systems made from two atoms in the 
same electronic state (or in some accidental situation) 
this coupling of COM and REL motion vanishes within 
the harmonic approximation. 

An even larger reduction of the energy is observed, if 
only the anharmonicity is considered as is reflected by 
£^(6) in which the coupling of COM and REL motion 
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TABLE I: The effect of the trapping potential on the en- 
ergy spectrum of the 1st trap-induced state of *^Rb-*°K for 
the trapping parameters of the experiment reported in [l^ . 
The A values defined by Eqs. p2|l - p4p are given in units of 
ft" ^ /kHz and calculated at asc = 6500 ao. 





^gcom 


^(2) 


^(6) 


Atot 


27.5 
40.0 


6.797 
6.828 


3.022 
3.689 


5.665 
6.243 


12.462 
13.071 



is ignored. For the considered system the effect of an- 
harmonicity is thus a larger correction to the separable 
harmonic approximation than the one caused by the cou- 
pHng of COM and REL motion. A further energy reduc- 
tion is found, if both effects are considered which leads to 
5^^^ . Interestingly, the energy reduction indicated by 
is larger than the sum of the energy reductions obtained 
separately for f and E''^^. The coupling of COM and 
REL motion is thus enhanced, if the more realistic sextic 
potential is considered instead of the harmonic one. 

In order to quantitatively describe the different effects 
of the trapping potential the energy differences 

Ageom = i?i'^-i?f^ (12) 

AiX = (13) 

Atot = i?f ' - = Ageom + A(6) p (14) 

may be introduced, where i = {lb, Iti}. Agoom charac- 
terizes the effect of anharmonicity of the optical lattice 
based on the uncoupled solutions. AcSup is a measure of 
the coupling between COM and REL motion within the 
harmonic {n ~ 2) or sextic (n ~ 6) potential. Finally, 
Atot specifies the energy difference between the simple 
harmonic approximation (in which the coupling of COM 
and REL motion is ignored) and the full solution of two 
atoms in a single site of an optical lattice (within the 
sextic approximation). 

As is evident from Fig. [21 the effect is largest for the 
strongly repulsive regime. This is due to a rise of the 
energy level to the region of higher anharmonicity of the 
trapping potential. Moreover, the state in this point is 
also long-ranged due to the strong repulsive interaction. 
Numerical values of the differences A ^21114^ for ^^Rb- 
^''K , the experimental trap parameters in [l3| . and = 
6500 ao are given in Table [D 

For the considered system the value Atot and thus 
the total energy between the uncoupled and the cou- 
pled harmonic approximation amounts to about 13 kHz 
for flsc — 6500ao. An effect of this size should be visi- 
ble in the experiment in [ist with a claimed resolution 
of 1.7 kHz but would not be resolvable with a ten times 
worse resolution as it occurs for a ten times shorter rf- 
pulse as was used, e.g., in [T^ . 

In an analysis of the influence of the interaction 
strength as is performed in this work it is important to 
stay within the restrictions of a single-site model. The 



parameter variation has to avoid situations in which tun- 
neling or even over-the-barrier transfer of atoms between 
different sites of a physical optical lattice can occur, since 
this range is clearly not adequately described with a har- 
monic or sextic potential with infinite walls. For exam- 
ple, for a very large positive scattering length the large re- 
pulsive interaction shifts the lowest-lying atom-pair state 
(1st trap-induced state) above the barrier of a true op- 
tical lattice. While this physical lattice would not sup- 
port any bound states, the harmonic or sextic potentials 
would still possess an infinite number of them. As is dis- 
cussed in Sec. lIII Cl it was always checked that the wave 
functions remain well localized within the boundaries of 
a single site of the optical lattice for the parameters used 
in this study. 



B. Other systems 

In order to obtain a more complete picture of the an- 
harmonicity and the coupling effects other systems may 
be analyzed. Besides the already considered ^^Rb-'^'^K 
pair (example of large masses and polarizabilities) other 
experimentally relevant alkali metal dimers like ^Li-^^'^Cs 
(small mass and polarizability of ^Li and for ^"^^Cs both 
characteristics are large) and ^Li-^Li (small masses and 
polarizabilities of both elements) arc investigated. 

Figure [3] shows the differences A (|12m4p as a function 
of the scaled interaction parameter ^ (see Sec. Ill D[) for 
different lattice depths obtained by the laser intensity 
variation for the three mentioned systems. As is evi- 
dent from Fig. [31 the harmonic coupling difference Acoup 
is not influenced by the lattice depth, because the cou- 
pling depends only on the polarizabilities and the masses 
and is of the form (/i2Q!i — I^ict2)- Therefore, Acoup is 
largest for ^Li-^^-^Cs and smallest for ^Li-^Li as is clear 
from Fig. [3] Beyond the harmonic approximation the 
mass, polarizability, laser intensity, and fc-vector depen- 
dence are mathematically non-trivial in the framework of 
the present approach. As a result, the different A val- 
ues have a behavior which is difficult to predict a priori. 
For example, while the total difference Atot decreases 
with the laser intensity for '^^Rb-^^K and increases for 
the other two systems, the values Ag and Agoom change 
their behavior not only with the laser intensity but also 
when going from one dimer to the other. Most notewor- 
thy, the A values for ^Li-''Li are not smaller than for the 
other pairs although this system is almost homonuclear. 

Another peculiar feature of the ^Li-^Li dimer com- 
pared to the other considered ones is the occurrence of 
negative values for Acoup in the case of large positive 
values of ^ and the laser intensity of 60Wcm^^. This 
leads to a smaller value of Atot compared to Agoom f^or 
these parameters. Clearly, the conclusions obtained for 
the generic '^^Rb-^^K system are not always transferable 
to other alkali metal dimers. 
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FIG. 3: (Color online) The energy differences [see Eqs. I|12I14^ ] 
Acolp (blue dots), A^o'up (red dashes), Agoom (green chain), Atot 
(black solid) (in multiples of both in atomic units) for different 
alkali metal dimers and intensities of the lattice laser (as specified 
in the graphs). The wavelength of the trap laser is 1030 nm. (The 
laser intensity 200 W/cm^ corresponds to a 30 lattice depth 
for the 87Rb-*'^K dimcr.) 



C. Wave- function analysis 
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FIG. 4: (Color online) Radial pair densities of the 1st trap-induced 
state of the *^Rb-''OK system in a 3D cubic lattice of 40 E^"" depth 
and the laser wavelength A of 1030 nm in the uncoupled harmonic 
(black solid), the harmonic with coupling (red dashes), the uncou- 
pled sextic (blue dots), and the sextic with coupling (green chain) 
approximations for the three interaction regimes: a) strongly at- 
tractive, b) almost zero interaction, c) strongly repulsive. (The 
insets show the densities for the least bound-state.) 



1. Radial pair densities 

An alternative analysis of the anharmonicity and 
COM-REL coupling effects is possible from the wave 
functions of the first trap induced and the least bound 
state. Since the probability density for finding a two- 
particle separation to lie in between r and r + dr is de- 



termined by the radial pair density 




it is convenient to discuss radial pair densities instead of 
the wave functions. In Eq. (jTS]) the function \xi) stands 
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for or |4>i) depending on the considered approxima- 
tion, dVfj is the COM volume element, and fir is the 
angular part of the REL motion coordinates. 

The energy spectrum of the first trap induced and the 
least bound states for the wide range of the interaction 
regimes was presented in Fig. [2] However, the three 
asymptotic interaction situations, namely, strong attrac- 
tion (csc — >■ ~oo), the almost zero interaction {age — > 0) 
and the strong repulsion (age +oo) arc foimd sufficient 
for the wave function analysis. Figure [4] shows the radial 
pair densities at the different levels of approximation for 
the three interaction regimes. As is evident from Fig.|4l a 
large attractive interaction leads to a very confined func- 
tion for the first trap-induced bound state while a large 
repulsive interaction does not only result in a node but 
also in a shift of the outermost lobe to larger interatomic 
distances. This shift is counteracted by the confinement 
of the trap. Remind, for the trap-free situation in the 
strongly repulsive regime the wave function crosses the 
intcrnuclcar axis exactly at the value of the scattering 
length. Evidently, the behavior of the density for almost 
zero interaction is only determined by the trap. 

As is apparent from Fig. [4] the inclusion of the an- 
harmonicity and the coupling leads to more extended 
pair densities. For all interaction regimes the following 
behavior is found. The effect is smallest for the har- 
monic coupling correction. A larger effect is found for 
the scxtic non-coupled case which is strengthened by the 
scxtic coupling. Effectively, the particles experience a 
more extended trap, if a more complete description of 
the problem is achieved. Such an effect is expected for 
the harmonic-to-sextic uncoupled description, since the 
harmonic trap is tighter than the sextic one as is evident 
from the sketch in Fig. [TJ While this is expected for the 
inclusion of the anharmonicity, this is not immediately 
clear for the coupling. 

The least bound state in the strongly repulsive regime 
is very long-ranged as the inset of Fig. \^c) shows. This 
distance is almost comparable with the extension of the 
first trap- induced state for the strongly attractive regime. 
Nevertheless, the influence of the anharmonicity and the 
coupling on the least bound state is almost absent, be- 
cause the state is energetically deeply bound and there- 
fore does almost not probe the anharmonicity of the lat- 
tice. 

In the considered parameter ranges the radial pair den- 
sities approach zero clearly before the interatomic dis- 
tance reaches the boundary of a single lattice site, i. e. for 
r < A/2, as can be seen from Fig. 2) Therefore, tunneling 
or a distribution of the dimer over more than a single lat- 
tice site does not occur and the present single-site model 
is applicable. Furthermore, effects of the artificial infi- 
nite walls of the harmonic or sextic potentials should not 
be a problem. However, the radial pair densities provide 
only an indication for the applicability of the single-site 
approximation, since it is still possible that the dimer 
as a whole may be distributed over more than one site. 
This can only be excluded form an analysis of the total 



wave function including the COM motion as is done in 
the next section. 

For the present type of the coupled problem analysis 
the plots of the radial pair densities can be used as a 
good check of the validity of the pseudopotential approx- 
imation. Varying the value of the scattering length in the 
pseudopotential approach it is possible to match the func- 
tion of the full solution at the long range distance. Such 
an approach is comparable to the energy-dependent con- 
cept developed for two identical particles in a harmonic 
trap [l^ and will be discussed in a separate work. Before 
entering such an analysis, it may be remarked that the 
radial pair densities as shown in Fig. [4] are also of interest 
for the investigation of the validity of a pseudopotential 
approximation for the interatomic interaction. In fact, 
it may even be used in order to attempt to obtain an 
improved pseudopotential describing atomic pairs in an 
optical lattice. A corresponding study is presently un- 
derway. 

2. Wave function in absolute coordinates 

It is instructive to analyze the full wave function or cor- 
responding particle density also in absolute coordinates 
of the laboratory space (ABS). They supply the complete 
information about the dimer and provide pictures of the 
COM and REL motion simultaneously. This is evidently 
not the case for the radial pair density (Fig. [?]) that is 
averaged over the COM motion and thus does not reveal 
whether the pair as a whole moves through the lattice. 
Note, while the radial pair density provides nevertheless a 
rather easy to interpret picture of the underlying physics, 
this is far less the case for its angular part. The reason 
is that the anisotropic (egg-box like) shape of the (cubic) 
optical lattice does not trivially show up in the REL coor- 
dinate system. An analysis in the lab frame is, however, 
also non-trivial, since the functions depend on six spa- 
tial coordinates. Use of the cubic symmetry reduces the 
size of the symmetry non-equivalent space, but it is still 
impractical to consider the complete multidimensional 
function. Instead, some insight may be gained from se- 
lected cuts. Although a number of cuts was analyzed in 
this work, only the results for cuts along the x coordinate 
of both atoms and thus for yi ~ Zi ~ (for both particles 
i) are shown and discussed. 

In order to quantitatively describe the different effects 
of the trapping potential it is again useful to consider 
not the wave functions at different level of approximation 
themsclf, but their respective differences. Similarly to 
the energy differences defined in Eqs. p2m4p . the wave- 
function differences 

AFgeom(ri,f2) = $r(fi,?^2)-$f^(ri,f2) (16) 

AFi")p(fi,f2) = 1>l")(rl,f2)-*f)(fi,f2) (17) 

Ai^tot(ri,r2) = $f^(fi,f2)-*f (fi,f2) (18) 

may be introduced. Cuts through these difference func- 
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FIG. 5; (Color online) Cuts of AF (defined in Eqs. H16ll8| l') along 
the X direction (yRb = Vk = -^Rb = zk = 0) for the 1st trap- 
induced state and almost non-interacting ®'^Rb-'"'K atoms in a 
3D cubic lattice (40 E^^ , A = 1030 nm): a) AFgcom(a;H,b, a^x) b) 

AFio^p(xRb,XK) c) AFc'oup(2:Rb,a;K) d) AFtot(a;Rb, i^k)- The dif- 
ferences A_F are given in atomic units, downscaled by the corre- 
sponding factors given in red. The black lines indicate the COM 
axis (xRb = — AtK/^Rb ^k)- The insets show \ AF\ (on an enlarged 
scale). 



tions AF are shown in Fig. [5] for the first trap-induced 
state and the almost non-interacting case (osc ~ 0). The 
sign convention used in Eqs. pGlllSp means that positive 
maxima in Fig. [5] correspond to the case that the wave 
functions in lower order of approximation have a larger 
amplitude than those in the higher one. (This choice is, 
of course, arbitrary and basically motivated by the fact 
that it leads to positive maxima in the center of the plots 
which is more suitable for optical reasons.) 

The diagonal = xk defines the REL coordinate 
axis. The wavefunctions and therefore also their differ- 
ences AF are strictly zero along the REL axis, since the 
molecular interaction potential rises exponentially to in- 
finity for r — > 0. Note, even for Osc = the atoms interact 
in the present approach, since the scattering length char- 
acterizes only the effective long-range interaction. In the 
case of the often adopted (5-type pseudopotential descrip- 



tion the interaction vanishes completely for age = and 
the wavefunction does not vanish around r = 0. Slightly 
away from the REL axis the wavefunction shows rapid 
oscillations due to the nodal structure that is again a 
consequence of the realistic interatomic interaction po- 
tential used in the present work. They are, however, not 
resolved in Fig [51 as these oscillations occur in a very 
small 10"'^ A/2) r range compared to the one dis- 
played. 

The COM axis is defined by = — MK/^Rb a;K- Since 
**^Rb-''°K is heteronuclear, the COM axis is rotated from 
the a;Rb = —xk diagonal and is thus for better readability 
explicitly indicated in the graphs. Another consequence 
of the heteronuclear character is the elliptical shape of 
the figures contours that would be circular in the case of 
a homonuclear system. 

Figure ini^a) characterizes the geometrical effect of the 
anharmonicity of an optical lattice. Effectively, the sex- 
tic trap is more extended than the harmonic one. This 
leads to the decrease of the density at the center of the 
potential and an increased probability at the potential 
edges. Therefore, the probability to find Rb and K atoms 
at a larger distance from the center of the optical lat- 
tice is higher for an anharmonic trap compared with a 
harmonic one. Note, this probability redistribution is 
not homogeneous. For example, for xk ~ 0.3 A/2 and 
x^ib ~ — O.lA/2 a pronounced minimum of the function 
as is evident from Fig.[5](a). This behavior 
in the ABS space is a direct consequence of the different 
COM and REL motion trapping depths. The COM of 
the system is more confined. Hence, the density shift in 
the direction of the COM axis is larger than for the REL 
one, as is better seen in the inset of Fig. EJa). In the 
inset one notices also that there exist small minima at 
the places where both atoms are close together (r « 0), 
but away from the center of the lattice site. In general, 
AFgeom is rather symmetric with respect to the COM 
and REL axes. 

Figure \5lh) shows the effect of the coupling of COM 
and REL motion within the harmonic approximation. 
For its understanding it is important to keep in mind 
that the definition of coupling between the different de- 
grees of freedom depends on the adopted coordinate sys- 
tem. In the present work it is defined by the Hamiltonian 
in Eq. ^ and thus the couphng of COM and REL coor- 
dinates. While this is a natural choice for discussions of, 
e. g., the radial pair density, its meaning is less transpar- 
ent for a discussion of wave functions in ABS coordinates 
of the two atoms. An evident example is the case of two 
truly non-interacting atoms in a harmonic trap. Even 
for a heteronuclear atom pair the problem separates in 
ABS coordinates, as was already mentioned in Sec. Ill Al 
However, treating this system in COM and REL coordi- 
nates the coupling term W in Eq. ([3]) and thus also the 

(2) 

difference AFcoiip is non-zero, but the latter refiects the 
non-separability due to the adopted coordinate system. 

A comparison of AFcoup shown in Fig. E^b) with the 
one obtained from the analytically known harmonic so- 
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lutions in cither ABS or REL and COM coordinates for 
truly non-interacting particles (see, e.g., [s^l) confirms 
that the structures in Fig. EJb) for the ®^Rb-^°K dimer 
with the long-range interaction being tuned to be al- 
most vanishing are similar. In contrast to the case of 
the geometry effect visible in Fig. [5][a) the maxima and 
minima in Fig. [5jb) have a more similar magnitude (the 
maxima being about 16 % larger in absolute value than 
the minima) and, clearly, they also originate from the 
different effects discussed above. Furthermore, they are 
located away from the center of the optical-lattice site. 
In fact, they are found, if one of the two atoms is lo- 
cated closely to the center and the other one is separated 
by about the most likely separation (about 0.2 A/ 2, see 
Fig. |4] (b)). The maxima (minima) are connected with 
the lighter K (heavier Rb) atom being close to the cen- 
ter. The anti-clockwisc rotation of the maxima and min- 
ima around the origin is due to the hcteronuclcar char- 
acter and reflects the coupling term in ABS coordinates, 
namely, A^Loup ^ e-T^K^jRb (where 7 is some constant). 
Also the different widths of the maxima and minima is a 
consequence of the heteronuclear character of ^^Rb-'^^K 
. While the off-centered minima of AFgeom are centered 
on the COM axis, the COM axis appear to separate the 

(2) 

minima and maxima of Ai^coup, although it does not de- 
fine a strict nodal plane. 

The sextic coupling effect presented in Fig. \Slc) is sim- 
ilar to the harmonic one in Fig. [5fb). However, the two 
maxima are now almost connected (if there were not the 
strict node on the REL axis) and form more a kind of 
plateau. The minima are less pronounced and as a conse- 
quence, the absolute values of the maxima are now about 
40% larger than the ones of the minima. 

Figure[5ljd) presents the complete effect of anharmonic- 
ity and coupling of the optical lattice. Compared to 
the previously discussed AF functions the shown A^tot 
is in shape most similar to Ai^ioup in Fig. [DJc). How- 
ever, the two maxima at the corners of the plateau ap- 
pear now to be merged with the central peak due to 
the new scale (and are basically only separated by the 
node along the REL axis). As a consequence, the den- 
sity of the exact sextic solution is reduced at the cen- 
ter of the lattice compared with the uncoupled harmonic 
approximation. In fact, as is evident from the equality 

AFtot = AFgeom + AFioup (compare Eq. ([14])), the merg- 
ing of the two maxima is simply a consequence of the 

superposition of the structures of AFioup and AFgeom- 
Note the correspondingly almost by a factor 2 larger am- 
plitude of AFtot compared to the other wavefunction dif- 
ferences. Since the minima of AFgcom and AFioup ap- 
pear at rather different places, their relative importance 
diminishes in comparison to the maxima. This leads to 
an about 66 % larger absolute value of the maxima com- 
pared to the minima in the case of AFtot. However, the 
additivity leads to an effective broadening of the minima 
of AFtot in direction of the COM axis compared to the 
mmrnia found for AFioup. 



To conclude the almost non-interacting case, the opti- 
cal lattice is in the coupled sextic description effectively 
more extended than in the uncoupled harmonic one. As 
a consequence of this anharmonicity the wave-function 
amplitude at the center of the lattice site decreases and 
is redistributed to the edges of the potential. As a conse- 
quence of the coupling, the decrease of the wave-function 
amplitude stretches further out along a diagonal in be- 
tween the COM and REL axes close to the axis defined 
by the Rb atom being located at the center of the lattice 
site (sRb = 0). On the other hand, the coupling leads 
also to minima (increased amplitude) along a diagonal 
between the COM and REL axes, but close to the xk = 
axis. As a consequence of the heteronuclear character of 
the ^^Rb-^^K dimer, the two diagonals are rotated with 
respect to the two corresponding x = axes. 




FIG. 6; (Color online) Cuts through the differences A_F as 
in Fig. [5l but a) AFgcom (aiRb, a:K) and c) AFtot (3:Rb, 3;k) for 
strongly attractively interacting particles (asc = —6600 ao); b) 
AFgcom (a;Rb, 2:k) and d) AFtot (x^b, i^k) for strongly repulsively 
interacting particles {asc = +6600 ao). 

Figure [HI characterizes the anharmonicity and coupling 
effect in the strongly interacting regimes. The most evi- 
dent difference between the (almost) non-interacting case 
and the strongly interacting situations is the pronounced 
squeezing of the central peak along the COM axis. This 
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is easily seen by comparing the geometry effect charac- 
terized by AFgeom in Figs. [D^a) and (b) with Fig. [5lja). 
Connected with this squeezing is an increase of the max- 
ima by a factor of almost 6 (strong repulsion) or more 
than 8 (strong attraction). The additional lobe occur- 
ring for strongly repulsive interaction at large distances 
leads to two further maxima on the COM axis (one for a 
positive and one for a negative value of xk) for A^Fgcom 
(Fig.IHb)), indicating a corresponding difference between 
the uncoupled harmonic and sextic solutions that occurs 
also at the outer lobes. The mimima on the COM axis are 
in this case shifted to larger distances from the REL axis. 
This is not the case for a strong attractive interaction, 
but there the amplitude of the minimum is even smaller 
than in the attractive case where it is already of less rela- 
tive importance compared to the central maxima than in 
the non-interacting case. On the other hand, the minima 
on the REL axis that had been very weak compared to 
the ones on the COM axis for the non-interacting case are 
in the strongly interacting cases much more pronounced, 
but also squeezed into a narrow regime close to the REL 
axis. 

As for the almost non- interacting case (Fig. EUd)), the 
total differences Ai^tot for the strongly interacting cases 
(Figs. IHtc) and (d)) differ from their Ai^geom counter- 
parts by the occurrence of two minima along diagonals 
between the COM and REL axes for xk ~ 0. While the 
also coupling-induced maxima for x^h ~ lead for the 
non-interacting case to a broad central peak, they appear 
in the strongly interacting case as shoulders. The reason 
is the massive squeezing of the central peak already dis- 
cussed for AFgoom- The additional maxima along the 
COM axis in the case of strong repulsion lead to a rather 
structured difference surface Ai^tot in this case. Another 
interesting effect visible from Fig. [DJd) is the enormous 
increase of the central maximum when comparing A_Ftot 
with AFgeom- For both the almost non-interacting and 
the strongly attractive case there is an approximate in- 
crease by a factor of 2, but in the strongly repulsive case 
there is a factor of more than 6. 

Compared to the analysis of the radial pair densities in 
Sec. nil C II it is evident that the absolute wave-function 
analysis reveals much more subtle details. In the case 
of radial pair densities there was the clear trend that 
improving the level of description leads to an increas- 
ing shift of probability from the maxima towards large 
separations. Similarly, the energies were uniformly low- 
ered. (Remind, however, that the energy analysis for 
"Li^Li showed that such a uniform trend is not found 
for all hctcronuclcar systems.) The cuts through the full 
wave functions show that the effects of coupling and an- 
harmonicity are not as trivial. Most importantly, they 
indicate that there is a lot of changes of the wave func- 
tions for short internuclear separations where, e.g., a 
pseudopotential approach is questionable. The relative 
importance of this regime of interatomic separations is, 
however, reduced, if an average over the angles is per- 
formed; simply because it scales with the radial part of 



the volume element, r^. This is also the reason why the 
energies are not very sensitive to this short-range regime 
and thus the pseudopotential approach may rather suc- 
cessfully predict also energy differences between different 
levels of approximation. 

The wave functions (not their differences) were also 
used in order to assure that the parameters chosen in this 
work allow a discussion in terms of a single site of an opti- 
cal lattice. Different cuts through the wave functions (in 
different directions relative to the optical lattice) never 
indicated a substantial wave function amplitude close to 
the boundaries of the single lattice site. 



D. Comparison to experiment 

A natural application of the present approach is to 
model the experimental results of C. Ospelkaus et al. p^ . 
In that experiment rf association was used to create 
molecules from fermionic and bosonic ^^Rb atoms in 
a 3D cubic optical lattice. The binding energy of the het- 
eronuclear molecules was measured as a function of the 
strength of an applied magnetic field. Figure [7] shows the 
experimental data for a lattice with depth VRb = 40 E^^ 
and wavelength A — 1030 nm. 

Note, the binding energies measured in the experiment 
are not the usual ones. In free space, real molecules (RM) 
close to the Feshbach resonance exist only on the repul- 
sive side of the resonance (csc > 0). The binding en- 
ergy measured in a trap-free situation is the one relative 
to the threshold energy of the continuum. In the pres- 
ence of an external optical lattice this continuum is dis- 
cretized, and there is instead a first trap-induced state. 
On the attractive side of the resonance (age < 0) the 
energy of this state is lowered relative to the field-free 
position. This leads to confinement-induced molecules 
(CIM) [l^ . In the experiment reported in [l^l , the bind- 
ing energy of the RM and CIM in a trap were measured. 
For Osc > the excitation energy of the repulsively in- 
teracting bound pair (RIP) [111 where repulsion between 
bosons and fermions shifts the two-particle ground state 
towards a higher energy was also measured. The cor- 
responding RM, CIM and RIP branches are denoted in 
Fig.[7]in which the experimental results of [ist are repro- 
duced for comparison. 

In order to compare the experimentally measured bind- 
ing energies with the theoretically calculated ones a 
proper mapping must be applied. Figure [S] outlines the 
procedure of how the binding energies were determined in 
the model. The scattering length Obg — 185ao is cho- 
sen as the i^-field-free background scattering length. The 
energy of the first trap-induced state obtained with the 
full sextic solution at Obg is chosen as energy zero and is 
marked explicitely in Fig. [SI A variation of the scattering 
length leads to energy shifts of the least bound and the 
first trap-induced states relative to this energy zero. The 
binding energy is a function of this shift as is indicated 
by arrows in Fig. [S] Specifically, the binding energy as a 
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FIG. 7: (Color online) The experimentally measured binding en- 
ergy (diamonds) of *^Rb-*'^K in an optical lattice (40-5^'', A of 
1030 nm) together with the theoretically calculated ones for the 
sextic potential and the energy-independent (black solid) or energy- 
dependent (red solid) scattering length and the Feshbach resonance 
parameters Bq = 546.8 G and AB = —3 G. The figure also shows 
the binding energy for the sextic trap and the energy-independent 
scattering length for the alternative value = 546.66 G (green 
solid). (All theoretical curves are full CI solutions.) 



FIG. 8: (Color online) Sketch of the procedure for obtaining bind- 
ing energies from the model. The 1st trap-induced (blue) and least 
bound (red) state energies obtained with a full CI calculation for 
the sextic potential (the same as in Fig.[2{b)) are shown. While the 
energy-offset of the first trap-induced level (blue arrows) relative 
to the energy zero (£-[^j' (a]-,g)) is responsible for the confinement 
phenomena, the energy offset of the least bound state (red arrow) 
is responsible for the pure molecular ones. 



function of the scattering length may be obtained from 
the present theoretical data with help of the relation 

E^,")(a3c;^)=Ei:)(abg)-E|")(a,e) (19) 

where i and n are, as before, i = {lb, Iti} and n ~ {2, 6}. 
Furthermore, E stands for the state energy at a given 
level of approximation E = {E,£}. This definition of a 
binding energy results in three different branches. While 
the first trap-induced state is responsible for the RIP and 
the CIM branch, the least bound state is responsible for 
the RM part. The corresponding branches are indicated 
in Fig. [7] and in the sketch of Fig. [51 

Experimentally, the binding energies were measured as 
functions of the magnetic field while theoretically calcu- 
lated energies are functions of the interaction strength 
represented by the scattering length (see Sec. Ill Dl for de- 
tails). To provide a B dependence of the theoretical data 
ttsc is mapped onto the magnetic field using a two-channel 
approximation [3l[ by the aid of 

asc{B) = (^1 - (20) 

where AB is the resonance width and Bq is the resonance 
position. Equation ([20]) gives in turn for the B field as a 
function of Ogc 

B(a,c) = AS f 1 - — ) +Bo . (21) 

V "bg/ 



The flgc values obtained from theory are inserted into 
Eq. ((2T|) to determine the B dependence of the energy. 

Figure [7] shows the binding energy obtained from the 
full sextic solution f for the experimental parameters 
of the trap and magnetic field Feshbach resonance pa- 
rameters AB = -3 G [13 and Ba = 546.8 G [36|. As is 
evident from Fig. [7] the model does not perfectly agree 
with the experiment. Some possible reasons of the dis- 
agreement are discussed in the following paragraph. 



1. Possible reasons of deviation between theory and 
experiment 

Equation ((20)) is derived for the lattice-free situation 
under the assumption that the collision between two 
atoms can be approximated by a two-channel scatter- 
ing model. In general, as was already mentioned in 
Sec. IIIDl the correct theoretical description requires a 
multi-channel scattering treatment which in the present 
case would have to incorporate also the optical lattice. 
Moreover, the present model uses an "artificial" vari- 
ation of the scattering length (see Sec. IIID[) . and the 
flsc values obtained from this variation are the ones of a 
single-channel approach. 

Even assuming the validity of Eq. ([^D)) there is an- 
other important factor influencing the comparison of the- 
ory and experiment. The values of the scattering length 
Osc and Obg of Eq. (|20p are determined in a lattice-free 
situation. In the presence of a trap these values must 
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be revised and adjusted to the trap parameters. It was 
shown in [l^ that the use of an energy-dependent scat- 
tering length gives ahnost correct energy levels for 
two harmonically trapped atoms. The evaluation of 
requires to solve the complete scattering problem and 
thus can only be obtained from the knowledge of 
the solution for the realistic atom-atom interaction po- 
tential. Eventually, the trap-free values of the scatter- 
ing length Osc and Obg must be substituted by appropri- 
ate energy-dependent scattering length values. How- 
ever, the problem is that the energy-dependent scattering 
length approach is so far developed only for the harmonic 
approximation, for s-wave collisions, and the uncoupled 
problem. An anharmonic, for example "sextic", energy- 
dependent scattering length concept as well as any other 
extensions of it do so far not exist to the authors' knowl- 
edge. 

In view of the absence of an beyond the uncou- 
pled harmonic approximation, the following procedure 
was adopted. The energy-dependent values of the scat- 
tering length are obtained using a solution for the pseu- 
dopotcntial energy and valid for a harmonic trap [37| 



1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 




545 545.5 546 546.5 547 547.5 548 
B(G) 

FIG. 9: (Color online) As Fig. [T] but using a harmonic poten- 
tial, an energy-independent scattering length, and the Feshbach- 
resonance positions Bq = 546.8 G (black solid) or Bq = 546.66 G 
(red solid). 
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where F is a gamma function. The energy of the REL 
motion obtained for the harmonic trap are imposed into 
Eq. p2|) and the values are obtained. The new values 
of the scattering length obtained with this manipulation 
are used for the mapping of the binding energies of £^^^ 
with the help of Eq. (PP]) . Figure [7] shows the result of 
this procedure. As is seen from the figure the shift of the 
spectral curve for the case of the energy-dependent scat- 
tering length along the _B-axis is not big but the curve 
is shifted along itself for the RIP branch and is tilted for 
the other ones in immediate proximity to the resonance. 
This may be seen as an indication that the energy depen- 
dence of the scattering length (properly included) does 
not have a too big effect, but the approximate implemen- 
tation is certainly not conclusive and thus cannot exclude 
a possible importance. 

Another important reason of the mismatch between 
theory and experiment could be an insufficient knowl- 
edge of the resonance parameters [l^ • It turns out to be 
sufficient to change the center of the Feshbach resonance 
to the value Bq = 546.66 G to match the experimental 
and the theoretical data. A variation of Bq of this size 
is well within the experimental uncertainty with which 
the resonance parameters are known [s^ . The result ob- 
tained with this modified value of Bq is also shown in 
Fig. [T] Remarkably, if both parameters Bq and AB are 
used together to fit the experimental curve it leads to 
a larger error than if only the parameter Bq is varied 
(see the discussion in the following subsection, especially 



Fig. [To)) . While a variation of AB and Bq shifts the the- 
oretical data along the magnetic-field axis, the variation 
of Obg leads in addition to a shift along the energy axis, 
since it changes the B-ficld-frce energy zero. 

Finally, one may address the question whether despite 
the number of uncertainties the effect of the anharmonic- 
ity and coupling is visible in the experiment [T3| . Figure[S] 
shows the binding energies obtained from the harmonic 
approximation. While the harmonic approximation pre- 
dicts the binding energy of the repulsively-interacting- 
pair part of the spectrum correctly, for other parts it 
results in a disagreement. A variation of the MFR pa- 
rameters does not lead to a simultaneous matching of all 
spectral branches. Therefore, it is possible to conclude 
that in the experiment [3l effects of anharmonicity and 
coupling (and thus deviations from a simple uncoupled 
harmonic model) were very likely detected. 



2. Comparison to a previous theoretical study 

The effects of anharmonicity and coupling of COM and 
REL motion in a single site of an optical lattice were 
also the subject of a recent theoretical study by Deuret- 
zbacher et al. p^ . The approach therein differs from 
the present one, since (i) it docs not use the full inter- 
atomic interaction potential but resorts to the pseudopo- 
tential approximation, (ii) a different partitioning of the 
Hamiltonian is adopted, and (iii) different basis functions 
(eigensolutions of the harmonic oscillator) were adopted. 

The two independently developed approaches pro- 
vide the possibility to further check whether theory has 
achieved a sufficient accuracy to investigate the small de- 
viations from the simple uncoupled harmonic approxima- 
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TABLE II: Influence of difi'erent levels of approximation on 
the energy of the 1st trap-induced state for three heteronu- 
clear systems. All results are obtained for asc = 6500 ao (or 
C(RbK) = 3.34, C(LiCs) = 3.76, ^(LiLi) = 3.24), lattice 
depths of Vi = V2 = 10i5r,roi where -Er.roi ~ k^/i'^fJ-), and 
a wavelength A — 1000 nm. 



atom pair 




A2 


f{6) 

"-•Iti "--Iti 


Atot 


*^Rb-''°K [present] 


3.79 


-0.12 


-0.29 


-0.41 


[18] 


3.74 


-0.12 


-0.27 


-0.39 


''Li-i^^'Cs [present] 


2.93 


-0.38 


-0.22 


-0.60 


[18j 


2.88 


-0.35 


-0.22 


-0.57 


^Li-^Li [present] 


3.93 


-0.01 


-0.30 


-0.31 


[18] 


3.92 


-0.01 


-0.29 


-0.30 



tion claimed to be found in the experiment in [13j . A 
consequence of difference (i) between the two approaches 
is furthermore the ability to inve stig ate the adequacy of 
the pseudopotential adopted in [1^. As a consequence 
of (iii) the approach in [18| can only be applied to very 
deep lattices and an extension to multiple-site lattices 
or even to shallow lattices is not straightforward. The 
reason is the rather strong spatial confinement of the 
harmonic-oscillator solutions. As a consequence, it needs 
an impractical large number of basis functions in order to 
cover an extended spatial regime. Since anharmonicity 
and coupling effects are different for shallower lattices 
as is discussed in Sec. IIIIBI th e tunneling effects may 
also play an important role j38l. [sot . Within the present 
approach calculations for multiple-wells and shallow lat- 
tices are straightforward and were already recently per- 
formed [TtI . 

The spatial compactness of the harmonic-oscillator 
eigenfunctions is on the other hand evident from the con- 
vergency study with respect to the Taylor expansion of 
the optical lattice performed in [l^ ■ As was discussed in 
Sec. IIII Al such a study is not senseful, since, e.g., even- 
order expansions lead to unphysical continua. Clearly, 
only a basis that does not explore the corresponding 
regime of the configuration space does not show any signs 
of these unphysical continua. 

Table |TT] shows a comparison of some energies and en- 
er gy d ifferences obtained with the numerical approach 
in [i^ and the present one for a large positive scattering 
length (osc = 6500 ao)- The results obtained with the 
two approaches do not differ very much in the case of all 
three considered alkali-metal dimers. The agreement of 
the energy differences is overall slightly better than the 
one of the absolute energies. The comparison seems to 
confirm the proper numerical implementation of both nu- 
merical approaches. Most importantly, it demonstrates 
that for the calculation of energy shifts as well as anhar- 
monic and coupling effects in a single site of an optical 
lattice the pseudopotential approach remains valid; at 
least to a very good approximation. 

It is presently not possible to attribute the remaining 
differences to the different atomic interaction potential 



or some remaining numerical uncertainty. Note, the dif- 
ferent interaction potential influences the results in two 
ways. First, the S-type pseudopotential does not prop- 
erly account for the short-range part of the interaction. 
Second, the mapping of the energy to a corresponding 
interaction strength is different in the two approaches. 
In the pseudopotential approach the scattering length is 
simply a parameter that enters the interaction potential, 
while it is extracted from the resulting wave function as 
was described in Sec. IIIDI As a consequence, there is a 
finite range in which age can be varied within the present 
approach. 

Comparable to the present finding (see Sec. IIII DT|) the 
binding-energy spectrum of *^Rb-^°K calculated in llSll 
does not agree very well with the experimental one in [l3| , 
if the previously experimentally determined Feshbach- 
resonance parameters (Bq = 546.8 G, AB = —3 G) 
are used. The authors in [l^ proposed that with the 
aid of the calculation it is in fact possible to improve 
on the MFR parameters. Such a fit (with the energy- 
independent scattering length) yielded the new resonance 
parameters Bq = 546.669 G and AB = -2.92 G [H]. 
This has to be contrasted with the present fit that yields 
the new resonance position Bq = 546.660 G but an 
unchanged width {AB = — 3 G), as was discussed in 
Sec. IIIDI Thus there is a similar (though slightly lar ger) 
trend for Bq, but disagreement with the results in [18| 
with respect to AB. 

In view of the different fit results, it is important to 
investigate in more detail their sensitivity on the fit pa- 
rameters. The quality of the fit depends on the agreement 
between the calculated binding energy 

(£(6) (5)) and the 

experimental one (ii^exp)- It is thus given by the relative 



E,^^{B) - £^^\B) 



£^oxp(-B) 



(23) 



Figure [10] shows 5{B) for three sets of MFR parame- 
ters: (i) Bo = 546.660 G and = -3 G (optimal 
fit parameters, this work), (ii) i?o = 546.669 G and 
AB = —2.92 G (optimal fit parameters found in [l^), 
and (in) Bq = 546.660 G and AB = -2.92 G (optimal fit 
parameter found in this work for Bq, but AB from p^). 
As is evident from Fig. [TOl any variation of either AB or 
Bq from their optimal values results in an increased error 
for all energy branches and all magnetic fields. Clearly, 
the fit shows a well defined minimum and thus there is no 
ambiguity in the fit parameters as it could occur, c. g., in 
the case of very shallow minima where the outcome of the 
fit may be determined by small numerical inaccuracies. 

Provided the fit fidelity in [l^ is comparable to the 
present one, i.e., a fit with the binding energies cal- 
culated in [l8[ using the optimal fit parameters of the 
present work would disagree with the experiment in a 
similarly pronounced fashion as shown in Fig. 1101 it is 
presently impossible to conclude whether theory has al- 
ready reached the level of accuracy that is required for an 
improved determination of MFR parameters. While both 
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FIG. 10: (Color online) The relative error defined in Eq. ||23|I for 
the alternative Feshbach parameters Bq = 546.660 G, AB = —3 G 
(squares), Bo = 546.669 G, AB = -2.92 G (circles), and Bo = 
546.66 G, AB = -2.92 G (triangles). 



fits appear to indicate a smaller value of Bo compared to 
the one previously extracted from experiment, the devi- 
ation between both fits is only about half as small as the 
improvement claimed in p^ . Clearly, such a result is 
from a statistical point of view inconclusive. In the case 
of the width AB the present finding agrees even fully to 
the previously determined value and thus disagrees with 



the result of the fit in [ll] . 

In order to obtain a more conclusive result it is vital 
to investigate whether the differences between the results 
in [l^ and the present ones are solely due to the use of the 
pseudopotential approximation or the more realistic in- 
teratomic interaction potential in the two studies. If this 
were the case, the fit results of the present study should 
be regarded as the more accurate ones. Furthermore, this 
would be an important example for the need to consider 
the interatomic interaction on a more accurate level than 
the one provided by the pseudopotential approximation. 
Since the implementation of the pseudopotential is due 
to the singular behavior of the 6 function non-trivial in 
the context of the present approach, such an investiga- 
tion has to be postponed to a separate work. Clearly, 
more experimental data (for different heteronuclear sys- 
tems) would also be very important for gaining a deeper 
insight and it is hoped that the present work stimulates 
such experimental activities. 

Finally, there are two further uncertainties in the de- 
termination of the MFR parameters from a fit like the 
one in |18j or in the present work. They are related to 
the way in which the mapping of the theoretical data 
onto the magnetic field is performed. As already men- 
tioned, this mapping is usually based on the assumption 
of validity of Eq. (PT|) and thus on the assumption that 
the B-field mapping of the multichannel MFR can be 
performed based solely on a scattering-length variation. 
Even in this case there is, however, the problem of the 
proper determination of the energy-dependent scattering 
length in an optical lattice which is so far unknown. The 
use of extracted from the harmonic uncoupled ener- 
gies for the mapping of the the full sextic energy results 
effectively in a shift of the energy-independent curve, as 
is seen in Fig. [T] However, both the energy-dependent 
and energy-independent Ogc discussed in this work uti- 
lize the same harmonic energy curve ignoring also the 
coupling to the COM motion. How the situation would 
change, if for a non-harmonic solution would be used, 
is difficult to predict, since the other curves in Fig. [2] 
not only differ in shape, but are also shifted relative to 
each other and contain the COM part. The overall good 
agreement of the theoretical binding energies (with fitted 
MFR parameters) to the experimental data is of course 
very suggestive that these uncertainties have a small in- 
fiuence, but this may be a pure coincidence. 

Both, the investigation of the appropriateness of the 
_B-field mapping as well as the question of the possibil- 
ity to define an energy-dependent scattering length be- 
yond the uncoupled harmonic approximation require a 
theoretical approach for the treatment of two atoms in 
an optical lattice as the one presented in this work and 
is presently pursued. It should, however, be emphasized 
that these uncertainties only affect the analysis of dimers 
close to a (magnetic) Feshbach resonance or, in general, if 
the proper description of the atom-pair requires a multi- 
channel treatment. The results of the previous sections 
arc valid independently of these uncertainties. Different 
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interaction regimes are experimentally accessible within 
the validity regime of a single-potential-curve treatment 
even for the same dimer by considering different isotopes 
or electronic states. The simplicity of experimental tun- 
ability as is found for magnetic Feshbach resonances is 
then of course lost. 



IV. CONCLUSION AND OUTLOOK 

An approach which allows for a full numerical descrip- 
tion of two ultracold atoms in 3D optical lattice is de- 
veloped. A detailed analysis of anharmonicity and cou- 
pling of center-of-mass and relative coordinates in terms 
of energy values and wave functions was performed for 
hetcronuclear dimers in a single site of an optical lat- 
tice. It is explained, why such a single site is optimally 
described by a sextic potential, if a finite Taylor expan- 
sion is used. The effects of deviations from the harmonic 
approximation and of the coupling were quantified and 
analyzed for different hetcronuclear systems, confinement 
strengths and interatomic interaction regimes. The influ- 
ence of the lattice is found to be always much stronger 
for the first trap-induced state than for the least bound 
state. As a consequence, binding energies are modified 
by the lattice mainly by the modification of the first trap- 
induced state. 

While the energy deviations from the harmonic un- 
coupled approximation is for all three considered generic 
hetcronuclear dimers largest for strong repulsive inter- 
action, the relative size (and even sign) of the energy 
change due to coupling or anharmonicity varies for the 
different dimers. The same is true for the the influence 
of the trap depth. A deeper lattice can lead to smaller 
or larger energy differences between the harmonic uncou- 
pled or the full coupled solution. While the analysis of 
the radial pair densities shows that the lattice mainly 
influences the maxima located at large interatomic sepa- 
rations, the analysis of cuts through the wave functions in 
absolute coordinates reveals non-negligible changes also 
at short interatomic distances. This may have important 
consequences for the validity of pseudopotential approx- 
imations. 

The results of the present theoretical approach are also 
compared to a recent experiment in which the binding 
energies of ^^Rb-''^°K have been measured as a function 
of an external magnetic field tuned close to a magnetic 
Feshbach resonance. The assumptions necessary for such 
a comparison are carefully discussed. It is found that 
very good agreement between experiment and theory can 



only be reached, if the previously experimentally deter- 
mined resonance parameters are modified. Since this 
needed modification is within the error bars with which 
the parameters had been determined before, this is not 
only reasonable, but may even indicate the possibility 
to more accurately determine the width and position of 
magnetic Feshbach resonances in ultracold atomic gases, 
as was proposed recently in a comparable theoretical 
study. However, the resonance parameters determined in 
the previous study based on the pseudopotential approx- 
imation differ from the ones found in the present work. 
If this deviation is due to the pseudopotential approxi- 
mation is difficult to judge at this moment. If this were 
the case, the then found breakdown of the pseudopoten- 
tial approximation would, of course, be a very interesting 
finding. A further investigation is therefore of great in- 
terest, and the present work stimulates hopefully also 
further experimental work in this direction. Since the 
infiucnce of anharmonicity and coupling becomes more 
important for less deep optical lattices and for excited 
trap levels corresponding experiments like the ones in 
j38l . [40I l4l| are expected to provide further tests of the 
approach presented in this work. 

Since the present approach was rather generally for- 
mulated and implemented, it allows immediately for fur- 
ther investigations that have partly been started or are 
planned for the future. This includes the consideration of 
highly anisotropic, asymmetric (disordered), or multiple- 
well lattice geometries. First results for triple- well poten- 
tials have, e. g., recently been used for the determination 
of Bose-Hubbard parameters and an inve stig ation of the 
validity of the Bose-Hubbard model itself [ij] . The study 
of an interesting physical phenomena like, e. g., the trap- 
induced resonances [43.l43j are also planned. Further ex- 
tensions of the approach should also allow to study the 
case of a pair of atoms or molecules interacting by non- 
centric, e.g., dipolar interactions. Finally, it is planned 
to extend the method for studies of the time-dependent 
dynamics of atomic pairs in time- varying lattices. 



Acknowledgments 

The authors are grateful to K. Sengstock, K. Bongs 
and F. Deuretzbacher for valuable discussion and to 
the Stifterverband fiir die Deutsche Wissenschaft, to 
the Deutsche Forschungsgemeinschaft (within Sonder- 
forschungsbereich SFB 450), and to the Fonds der 
Chemischen Industrie for financial support. 



[1] M. H. Anderson, J. R. Ensher, M. R. Matthews, C. E. 

Wieman, and E. A. Cornell, Science 269, 198 (1995). 
[2] K. B. Davis, M. O. Mewes, M. R. Andrews, N. J. van 

Druten, D. S. Durfee, D. M. Kurn, and W. Ketterle, 

Phys. Rev. Lett. 75, 3969 (1995). 



[3] D. Jaksch, C. Bruder, J. I. Cirac, C. W. Gardiner, and 
R ZoUer, Phys. Rev. Lett. 81, 3108 (1998). 

[4] M. Greiner, O. Mandel, T. Esslinger, T. Hansch, and 
I. Bloch, Nature 415, 39 (2002). 

[5] M. Kohl, H. Moritz, T. Stoferle, K. Giinter, and 



18 



T. Esslinger, Phys. Rev. Lett. 94, 080403 (2005). 
[6] D. R. Meacher, Contemporary Physics 39, 329 (1998). 
[7] I. Bloch, Nature Physics 1, 23 (2005). 
[8] M. Lewenstein, A. Sanpera, V. Ahufinger, B. Damski, 
A. Sen, and U. Sen, Advances in Physics 56, 243 (2007). 
[9] T. Loftus, C. A. Regal, C. Ticknor, J. L. Bohn, and D. S. 
Jin, Phys. Rev. Lett. 88, 173201 (2002). 
[10] C. A. Regal, C. Ticknor, J. L. Bohn, and D. S. Jin, Na- 
ture 424, 3 (2003). 
[11] J. Hubbard, Proc. R. Soc. London A 276, 238 (1963). 
[12] T. Stoferle, H. Moritz, K. Giinter, M. Kohl, and 

T. Esslinger, Phys. Rev. Lett. 96, 030401 (2006). 
[13] C. Ospelkaus, S. Ospelkaus, L. Humbert, P. Ernst, 
K. Sengstock, and K. Bongs, Phys. Rev. Lett. 97, 120402 

(2006) . 

[14] T. Busch, B.-G. Englert, K. Rzazewski, and M. Wilkens, 

Found, of Phys. 28, 549 (1998). 
[15] E. L. Bolda, E. Tiesinga, and P. S. Julienne, Phys. Rev. 

A 71, 033404 (2005). 
[16] S. Grishkevich and A. Saenz, Phys. Rev. A 76, 022704 

(2007) . 

[17] P. Schneider, S. Grishkevich, and A. Saenz, 
arXiv:0903.3185 (2009), (submitted to Phys. Rev. 
A). 

[18] F. Deuretzbacher, K. Plassmeier, D. Pfannkuche, 
F. Werner, C. Ospelkaus, S. Ospelkaus, K. Sengstock, 
and K. Bongs, Phys. Rev. A 77, 032726 (2008). 

[19] D. Blume and C. H. Greene, Phys. Rev. A 65, 043613 

(2002) . 

[20] Z. Idziaszek and T. Calarco, Phys. Rev. A 71, 050701(R) 
(2005). 

[21] I. S. Lim, M. Pernpointner, M. Seth, J. K. Laerdahl, 

P. Schwerdtfeger, P. Neogrady, and M. Urban, Phys. Rev. 

A 60, 2822 (1999). 
[22] S. Rousseau, A. R. AUouche, and M. Aubert-Frecon, 

J.Mol. Spec. 203, 235 (2000). 
[23] W. T. Zemke, R. Cote, and W. C. Stwalley, Phys. Rev. 

A 71, 062706 (2005). 
[24] A. Derevianko, J. F. Babb, and A. Dalgarno, Phys. Rev. 

A 63, 052704 (2001). 
[25] S. G. Porsev and A. Derevianko, J. Chem. Phys. 119, 844 

(2003) . 



[26] F. Ferlaino, C. D'Errico, G. Roati, M. Zaccanti, M. In- 

guscio, G. Modugno, and A. Simoni, Phys. Rev. A 73, 

040702 (R) (2006). 
[27] P. Staanum, A. Pashov, H. Knockel, and E. Tiemann, 

Phys. Rev. A 75, 042513 (2007). 
[28] J. Goldwin, S. Inouye, M. L. Olsen, B. Newman, B. D. 

DePaola, and D. S. Jin, Phys. Rev. A 70, 021601(R) 

(2004). 

[29] S. Inouye, J. Goldwin, M. L. Olsen, C. Ticknor, J. L. 

Bohn, and D. S. Jin, Phys. Rev. Lett. 93, 183201 (2004). 
[30] M. Zaccanti, C. D'Errico, F. Ferlaino, G. Roati, M. In- 

guscio, and G. Modugno, Phys. Rev. A 74, 041605(R) 

(2006) . 

[31] T. Kohler, K. Coral, and P. S. Julienne, Rev. Mod. Phys. 

78, 1311 (2006). 
[32] C. A. Regal and D. S. Jin, Phys. Rev. Lett. 90, 230404 

(2003). 

[33] C. F. Fischer and M. Idrees, J. Phys. B: At.Mol.Phys. 
23, 679 (1990). 

[34] A. S. Dawydow, Quantenmechanik (VEB Deutscher Ver- 

lag der Wissenschaften, Berlin, 1967). 
[35] K. Winkler, G. Thalhammer, F. Lang, R. Grimm, 

J. Hecker, Denschlag, A. J. Daley, A. Kantian, H. P. 

Biichler, and P. ZoUer, Nature 441, 853 (2006). 
[36] S. Ospelkaus, C. Ospelkaus, L. Humbert, K. Sengstock, 

and K. Bongs, Phys. Rev. Lett. 97, 120403 (2006). 
[37] B. Deb and L. You, Phys. Rev. A 68, 033408 (2003). 
[38] S. Foiling, S. Trotzky, P. Cheinet, M. Feld, R. Saers, 

A. Widera, T. Miiller, and I. Bloch, Nature 448, 1029 

(2007) . 

[39] S. Trotzky, P. Cheinet, S. Foiling, M. Feld, U. Schnor- 
rberger, A. M. Rey, A. Polkovnikov, E. A. Demler, 
M. Lukin, and I. Bloch, Science 319, 5861 (2008). 

[40] T. Miiller, S. Foiling, A. Widera, and I. Bloch, 
Phys. Rev. Lett. 99, 200405 (2007). 

[41] L. Fallani, J. E. Lye, V. Guarrera, C. Fort, and M. In- 
guscio, Phys. Rev. Lett. 98, 130404 (2007). 

[42] M. Olshanii, Phys. Rev. Lett. 81, 938 (1998). 

[43] R. Stock, I. H. Deutsch, and E. L. Bolda, Phys. Rev. Lett. 
91, 183201 (2003). 



